function [index_result,error_time,mean_pz,data_splined] = window_adaptive_size(index_points_error,data)
% 窗口自适应改变大小，并返回相应结果
%   窗口的大小可以根据异常轨迹段的长短进行自适应改变
%此处的轨迹段长短是值时间，即轨迹段的时间跨度
%窗口的大小需要与该轨迹段的时间跨度相适应
%index_points_error 异常轨迹段对应的轨迹点的索引
%data  原始轨迹数据
%error_time  位置异常数据对应的时间值
%mean_pz 原始、剔除、插值等不同阶段的轨迹偏转角的平均值
%data_splined 当前轨迹异常点对应的插值数据
%index_result 根据范围确定的删除行号

% 窗口大小是基于异常轨迹段的长短，即异常点数量
%window_size=length(points_error_index);%窗口大小，其值直接决定参考轨迹段的长短

%% 自适应窗口大小，分前窗口和后窗口
%窗口大小根据异常轨迹段的时间差值确定
data_error =data(index_points_error,:);
time_reference=data_error.time(end)-data_error.time(1);%异常轨迹段的时间跨度值
index_points_error(1)

%查找特定的数据点
% if ismember(439,index_points_error)
% test=1;
% end


%异常轨迹段之前的窗口成为前窗口
time_bf_reference=data.time(index_points_error(1)-1); %异常轨迹段之前轨迹点的时间值作为基准时间
time_diff_bf_reference= time_bf_reference -data.time(1:index_points_error(1)-1);%异常轨迹点之前所有点与基准时间之间的时间差值
% 找到差值小于或等于 time_reference 的最小值的数量，即前窗口大小
window_size_bf = sum(time_diff_bf_reference <= time_reference);

%后窗口大小
time_af_reference = data.time(index_points_error(end)+1);
time_diff_af_reference=data.time(index_points_error(end)+1:end) -time_af_reference;
window_size_af = sum(time_diff_af_reference <= time_reference);

%% 根据窗口大小，计算前后窗口轨迹段，分经纬度、角度值
%应该是综合窗口中基于经纬度判断要删除的点、基于角度值判断要删除的点，两者重合的部分才是真正应该删除的点。
%前面轨迹段
index_bf_start=index_points_error(1)-window_size_bf;%前窗口轨迹段开始位置
%防止越界
if index_bf_start<=0
    index_bf_start=1;
end

index_bf_end =index_points_error(1)-1;%前窗口轨迹段结束位置
%防止越界
if index_bf_end>= size(data,1)
    index_bf_end=size(data,1);
end
index_bf_window = index_bf_start:index_bf_end;%前窗口轨迹点索引
data_bf_window = data(index_bf_window,:);%前窗口轨迹段数据
lonlat_bf_window=[data_bf_window.Lon,data_bf_window.Lat];%前窗口轨迹段经纬度

%后面轨迹段
index_af_start=index_points_error(end)+1;
if index_af_start >= size(data,1)
    index_af_start=size(data,1);
end

index_af_end=index_points_error(end)+window_size_af;
if index_af_end >= size(data,1)
    index_af_end=size(data,1);
end
index_af_window =index_af_start:index_af_end;%后窗口数据索引
data_af_window = data(index_af_window,:);
lonlat_af_window=[data_af_window.Lon,data_af_window.Lat];%后窗口轨迹段经纬度

%% 经纬度范围判定
lonlat_mean_bf=[mean(lonlat_bf_window(:,1)),mean(lonlat_bf_window(:,2))];%前窗口经纬度的均值
lonlat_mean_af=[mean(lonlat_af_window(:,1)),mean(lonlat_af_window(:,2))];%后窗口经纬度的均值

%% 使用原始的矩形区域时，无法处理可能出现的抛物线顶点问题，因此采用圆形范围区域
%前后窗口经纬度的均值可以左右范围矩阵的两个顶点
% polygon_lon=[lonlat_mean_bf(1) lonlat_mean_bf(1) lonlat_mean_af(1) lonlat_mean_af(1) lonlat_mean_bf(1) ];
% polygon_lat=[lonlat_mean_bf(2) lonlat_mean_af(2) lonlat_mean_af(2) lonlat_mean_bf(2) lonlat_mean_bf(2)];
%
% % index_delete={};%保存要删除的位置索引
% for j=1:length(index_points_error)
%     %0表示不在内部，1表示在内部,因此只要统计in和on中0的位置和数量，就可知道需要剔除的点
%     [isIn, isOn]=inpolygon(data.Lon(index_points_error), data.Lat(index_points_error),polygon_lon,polygon_lat);
% end
%
% %获得基于经纬度判定的异常数据的索引
% index_error_lonlat_org=index_points_error((~isIn) & (~isOn));

%% 圆形区域法进行判定
%绘制圆,注意坐标转换
% 计算圆心（中点）
[x_mean_bf,y_mean_bf] =latlon2Cartesian_point_set(lonlat_mean_bf);
[x_mean_af,y_mean_af]=latlon2Cartesian_point_set(lonlat_mean_af);
[data_error.Lon,data_error.Lat] =latlon2Cartesian_point_set([data_error.Lon,data_error.Lat]);


center =[(x_mean_bf+x_mean_af),(y_mean_bf+y_mean_af)]/2;% (lonlat_mean_bf + lonlat_mean_af) / 2;
% 计算圆的半径（直径的一半）
radius = norm([x_mean_bf,y_mean_bf] - [x_mean_af,y_mean_af]) / 2;
% 计算检查点到圆心的距离
distances2center = sqrt((data_error.Lon - center(1)).^2 + (data_error.Lat - center(2)).^2);%求解异常点集中每个点到圆心的距离
% 判断每个点是否在圆内或圆上
in_circle = distances2center <= radius;

%如果不在圆上或者圆内，则就是异常点
index_error_lonlat = index_points_error(in_circle==0);
%如果检测到首尾节点是错误的，那需要将其进行剔除
index_error_lonlat(index_error_lonlat==index_points_error(1))=[];
index_error_lonlat(index_error_lonlat==index_points_error(end))=[];

if isempty(index_error_lonlat)
    index_error_lonlat=NaN;
end

%% 航向变化率计算
%前窗口航向变化率
diff_time_bf_window = diff(data_bf_window.time);
diff_trueheading_bf_window = diff(data_bf_window.Course);
trueheading_bf_deleta = diff_trueheading_bf_window ./diff_time_bf_window;
trueheading_bf_mean=mean(trueheading_bf_deleta);

%后窗口航向变化率
diff_time_af_window = diff(data_af_window.time);
diff_trueheading_af_window =diff(data_af_window.Course);
trueheading_af_deleta =  diff_trueheading_af_window ./diff_time_af_window;
trueheading_af_mean = mean(trueheading_af_deleta);

%如果是按照最开始的思路，计算的是异常数据之间的变化率，当存在错误点时，改点与之前和之后点之间的变化率均应该存在错误
%这样无形中增加了判定的难度
%因此，应该是异常数据分别与前、后窗口中的相近点（前窗口最后的点、后窗口开始的点）之间的变化率来判定某点是否属于离群点，即错误点
%错误数据与前窗口末尾数据之间的时间差值和角度差值，称之为"前窗口边界差值"
time_diff_bf_error=data_error.time-data_bf_window.time(end);%前窗口末尾值与错误数据的时间差
trueheading_diff_bf_error=data_error.Course-data_bf_window.Course(end);
trueheading_bf_error_deleta= trueheading_diff_bf_error./time_diff_bf_error;%前窗口边界航向变化率

%错误数据与后窗口开端数据之间的时间差值和角度差值，称之为"后窗口边界差值"
time_diff_af_error=data_af_window.time(1)-data_error.time;
trueheading_diff_af_error=data_af_window.Course(1)-data_error.Course;
trueheading_af_error_deleta= trueheading_diff_af_error./time_diff_af_error;%后窗口边界航向变化率


% %写文章可以作为示意图使用
% if ismember(54,index_points_error)
% figure
% yline(trueheading_bf_mean, 'b-', 'LineWidth', 2); % 前窗口航向平均变化率
% hold on
% yline(trueheading_af_mean, 'm--', 'LineWidth', 2);  % 后窗口航向平均变化率
% hold on
% plot(trueheading_bf_error_deleta,'k-*', 'LineWidth', 1);
% hold on
% plot(trueheading_af_error_deleta,'k-s', 'LineWidth', 1);
%
% legend('前窗口航向平均变化率','后窗口航向平均变化率','前窗口边界航向变化率','后窗口边界航向变化率');
%
% xticks(1:4) % 设置横坐标的刻度位置
% xticklabels({'P_1', 'P_2', 'P_3', 'P_4'}) % 设置横坐标的刻度标签
% xlabel('轨迹点');
% ylabel('航向变化率(°/s)');
% % 显示网格
% grid on;
% end




%% 航向变化率判定
%基于上述前窗口航向变化率均值和后窗口航向变化率均值，分别得到并绘制两条线
%根据窗口边界航向变化率值与对应的均值之间的差值的和来判断哪个点异常
%注意：异常点数据由于前期组合时存在扩张，因此只需要判断中间两个点是否异常即可
distance_bf = abs(trueheading_bf_error_deleta(2:end-1)-trueheading_bf_mean);
distance_af = abs(trueheading_af_error_deleta(2:end-1)-trueheading_af_mean);

distance_error = distance_bf + distance_af;
[~,max_pos]=max(distance_error);%最大值所在位置

index_error_trueheading =index_points_error(max_pos+1);%需要+1转换为原始位置
if isempty(index_error_trueheading)
    index_error_trueheading=NaN;
end

%% 以经纬度范围判定的数值为基准
%根据测试，经纬度范围判定的比角度值判定的准，因此逻辑是以经纬度范围确定数值后；
%上述经纬度和航向变化率的判定结果存在几种可能：
% 1、两者的判断结果是一致的，且都是实际数值
% 2、两者其中一个是空集，另外一个是实际数值
% 3、两者都是实际数值，但是数值冲突
% 第1和2处理方法：两者求合集并删除空集即可
%第3的处理方法：两者冲突时，以经纬度范围判定的结果为准


%经纬度判定结果和航向变化率判定结果求合集
index_result =union(index_error_lonlat, index_error_trueheading);

if any(isnan(index_result))%其中含有空集
    % 使用 isnan 函数生成逻辑数组，指示哪些元素是非空的
    nonEmpty = ~isnan(index_result);
    % 使用逻辑索引删除空数组
    index_result = index_result(nonEmpty);
else %没有空集
    if length(index_result)>1%第1\2情形
        index_result = index_error_lonlat;
    end
end

%判断错误类型，如果是位置异常数据（时间不重复），则需要插值
%则记录该错误数据的time值即可，以便后续在原始数据中查找对应的行
%如果此时错误数据中含有与该异常值相同的时间，则说明是部分重复数据类型，不需要差值，不记录，否则就记录
index_result_len=length(index_result);
error_time =zeros(1,index_result_len);%初始值为0，因此，如果时间time值为0的话，就是时间重复类型不需要差值
for i=1:index_result_len
    flag_error=abs(data_error.time-data.time(index_result(i)))<2; %如果时间差是0或者1，则认定为时间重复，由于前面进行了重复只剔除，因此差值为1的情形原则上被剔除了
    if sum(flag_error==1)==1  %如果没有时间重复，则是位置异常
        error_time(i) =data.time(index_result(i));
    end
end


%对应位置进行插值
index_result_len=length(index_result);
%进行相应的插值,插值计算必须放在此处的原因是，插值是基于窗口数据的
%只是对位置异常点进行插值，如果是部分重复数据，即如果时间相同，但是经纬度不同，则这种类型不插值
time_0='00:00:00';%时间完全相同

% mean_pz_struct =struct('org',-1,'deleted',-1,'spline',-1);%分别保存原始、剔除、插值等不同阶段的平均偏转角
% %但是需要注意的是，当index_result>1时，默认值需要对应着数量进行更改
% mean_pz=repmat(mean_pz_struct, 1, index_result_len);
mean_pz=zeros(index_result_len,3)-1;%分别保存原始、剔除、插值等不同阶段的平均偏转角
%错误类型,为后续是否进行插值进行提前预判
flag_error =  strcmp(char(data_error.Time(3)-data_error.Time(2)),time_0);%如果是位置异常数据，则是0类型，如果是部分时间重复类型，则是1类型。
data_error =data(index_points_error,:);%前面的异常数据的经纬度已经转换，因此需要重新赋值
data_window=[data_bf_window;data_error;data_af_window];
data_splined=[];
if ~flag_error
    data_splined = data(index_result,:);%异常点原始数据，需要用插值数据替换，

    data_window_deleted=data(setdiff(index_points_error,index_result),:);%异常数据集剔除异常行之后的数据
    data_window_deleted=[data_bf_window;data_window_deleted;data_af_window];%异常数据集（剔除异常点）+前、后窗口数据,注意顺序
    index_interp = data.time(index_result)-data_window_deleted.time(1);%需要插值点与数据起始点之间的时间插值,异常点可能不只一个，因此插值也不只一个
    for i=1:length(index_interp)
        [data_splined.Lon(i),data_splined.Lat(i),data_splined.Course(i)]= spline_process(data_window_deleted,index_interp(i));

        %分别计算对应的平均偏转角
        pz_org=calculate_track_degree([data_window.Lon,data_window.Lat]);
        mean_pz(i,1)=mean(pz_org);

        pz_deleted=calculate_track_degree([data_window_deleted.Lon,data_window_deleted.Lat]);
        mean_pz(i,2)=mean(pz_deleted);

        data_window_spline=[data_window_deleted;data_splined];
        data_window_spline=sortrows(data_window_spline,'time');
        pz_spline=calculate_track_degree([data_window_spline.Lon,data_window_spline.Lat]);%插值
        mean_pz(i,3)=mean(pz_spline);

    end

end

end


